Extended Data Figure 8d-f - Enrichment analysis
# Find markers of main populations
gc()
## used (Mb) gc trigger (Mb) max used (Mb)
## Ncells 11460755 612.1 20067912 1071.8 20067912 1071.8
## Vcells 3620957988 27625.8 6773128668 51674.9 6521432058 49754.6
options(future.globals.maxSize= 8*1024^3)
s.markers <- DElegate::FindAllMarkers2(macro, replicate_column = "patient", group_column = "anno_l1")
# universe
background <- rownames(macro)
ego_h <- list()
# pathways sets
GO_C5_gene_sets = msigdbr(species = "human", category = "C5", subcategory = "BP")
msigdbr_t2g_C5 = GO_C5_gene_sets %>% dplyr::distinct(gs_name, gene_symbol) %>% as.data.frame()
GO_H_gene_sets = msigdbr(species = "human", category = "H")
msigdbr_t2g_H = GO_H_gene_sets %>% dplyr::distinct(gs_name, gene_symbol) %>% as.data.frame()
CBioPlanet_gene_sets = read.csv2(params$path_to_CBioplanet, sep = ",")
CBioPlanet_gene_sets <- CBioPlanet_gene_sets[,c(2,4)]
colnames(CBioPlanet_gene_sets) <- c("gs_name", "gene_symbol")
CBioPlanet = CBioPlanet_gene_sets %>% dplyr::distinct(gs_name, gene_symbol) %>% as.data.frame()
# FNC1
signature_h <- s.markers$feature[s.markers$group1 == "FCN1" & s.markers$padj < 0.05 & abs(s.markers$log_fc)>0.58]
##C5
ego_h[["FCN1"]] <- enricher(gene = signature_h, universe = background, TERM2GENE = msigdbr_t2g_C5)
b <- barplot(ego_h[["FCN1"]], title = "MSigDB C5 enrichments in FCN1", showCategory = 4 )
tmp1 <- b$data
tmp1$ID <- str_trunc(tmp1$ID, 150, "right")
tmp1$ID <- factor(tmp1$ID, levels = tmp1$ID[order(tmp1$Count)])
tmp1$type <- "FCN1"
##Hallmarks
ego_h[["FCN1"]] <- enricher(gene = signature_h, universe = background, TERM2GENE = msigdbr_t2g_H)
b <- barplot(ego_h[["FCN1"]], title = "MSigDB H enrichments in FCN1", showCategory = 3 )
tmp2 <- b$data
tmp2$ID <- str_trunc(tmp2$ID,150, "right")
tmp2$ID <- factor(tmp2$ID, levels = tmp2$ID[order(tmp2$Count)])
tmp2$type <- "FCN1"
##CBioPlanet
ego_h[["FCN1"]] <- enricher(gene = signature_h, universe = background, TERM2GENE = CBioPlanet)
b <- barplot(ego_h[["FCN1"]], title = "CBioPlanet enrichments in FCN1", showCategory = 4 )
tmp3 <- b$data
tmp3$ID <- str_trunc(tmp3$ID, 150, "right")
tmp3$ID <- factor(tmp3$ID, levels = tmp3$ID[order(tmp3$Count)])
p1 <- ggplot(tmp1, aes(Count, ID)) +
geom_bar(stat = "identity", color="black", fill = colors_anno["FCN1"]) +
theme_classic()+
geom_text(
aes(label = (paste( "padj=", format(p.adjust,scientific = TRUE, digits = 2)))),
color = "black",
size = 6,
hjust=1,
position = position_dodge(0.5) ) +
theme(text = element_text(size=20)) +
ggtitle("FCN1") +
ylab("GO:BP")
p2 <- ggplot(tmp2, aes(Count, ID)) +
geom_bar(stat = "identity", color="black", fill = colors_anno["FCN1"]) +
theme_classic()+
geom_text(
aes(label = (paste( "padj=", format(p.adjust,scientific = TRUE, digits = 2)))),
color = "black",
size = 6,
hjust=1,
position = position_dodge(0.5) ) +
theme(text = element_text(size=18)) +
ggtitle("") +
ylab("Hallmarks")
# C1QC
signature_h <- s.markers$feature[s.markers$group1 == "C1QC" & s.markers$padj < 0.05 & abs(s.markers$log_fc)>0.58]
##C5
ego_h[["C1QC"]] <- enricher(gene = signature_h, universe = background, TERM2GENE = msigdbr_t2g_C5)
b <- barplot(ego_h[["C1QC"]], title = "MSigDB C5 enrichments in C1QC", showCategory = 4 )
tmp4 <- b$data
tmp4$ID <- str_trunc(tmp4$ID, 150, "right")
tmp4$ID <- factor(tmp4$ID, levels = tmp4$ID[order(tmp4$Count)])
tmp4$type <- "C1QC"
##Hallmarks
ego_h[["C1QC"]] <- enricher(gene = signature_h, universe = background, TERM2GENE = msigdbr_t2g_H)
b <- barplot(ego_h[["C1QC"]], title = "MSigDB H enrichments in C1QC", showCategory = 3 )
tmp5 <- b$data
tmp5$ID <- str_trunc(tmp5$ID,150, "right")
tmp5$ID <- factor(tmp5$ID, levels = tmp5$ID[order(tmp5$Count)])
tmp5$type <- "C1QC"
##CBioPlanet
ego_h[["C1QC"]] <- enricher(gene = signature_h, universe = background, TERM2GENE = CBioPlanet)
b <- barplot(ego_h[["C1QC"]], title = "CBioPlanet enrichments in C1QC", showCategory = 4 )
tmp6 <- b$data
tmp6$ID <- str_trunc(tmp6$ID, 150, "right")
tmp6$ID <- factor(tmp6$ID, levels = tmp6$ID[order(tmp6$Count)])
tmp6$type <- "C1QC"
p4 <- ggplot(tmp4, aes(Count, ID)) +
geom_bar(stat = "identity", color="black", fill = colors_anno["C1QC"]) +
theme_classic()+
geom_text(
aes(label = (paste( "padj=", format(p.adjust,scientific = TRUE, digits = 2)))),
color = "black",
size = 6,
hjust=1,
position = position_dodge(0.5) ) +
theme(text = element_text(size=18)) +
ggtitle("C1QC") +
ylab("")
p5 <- ggplot(tmp5, aes(Count, ID)) +
geom_bar(stat = "identity", color="black", fill = colors_anno["C1QC"]) +
theme_classic()+
geom_text(
aes(label = (paste( "padj=", format(p.adjust,scientific = TRUE, digits = 2)))),
color = "black",
size = 6,
hjust=1,
position = position_dodge(0.5) ) +
theme(text = element_text(size=18)) +
ggtitle("") +
ylab("")
p6 <- ggplot(tmp6, aes(Count, ID)) +
geom_bar(stat = "identity", color="black", fill = colors_anno["C1QC"]) +
theme_classic()+
geom_text(
aes(label = (paste( "padj=", format(p.adjust,scientific = TRUE, digits = 2)))),
color = "black",
size = 6,
hjust=1,
position = position_dodge(0.5) ) +
theme(text = element_text(size=18)) +
ggtitle("") +
ylab("")
# SPP1
signature_h <- s.markers$feature[s.markers$group1 == "SPP1" & s.markers$padj < 0.05 & abs(s.markers$log_fc)>0.58]
##C5
ego_h[["SPP1"]] <- enricher(gene = signature_h, universe = background, TERM2GENE = msigdbr_t2g_C5)
b <- barplot(ego_h[["SPP1"]], title = "MSigDB C5 enrichments in SPP1", showCategory = 4 )
tmp7 <- b$data
tmp7$ID <- str_trunc(tmp7$ID, 150, "right")
tmp7$ID <- factor(tmp7$ID, levels = tmp7$ID[order(tmp7$Count)])
tmp7$type <- "SPP1"
##Hallmarks
ego_h[["SPP1"]] <- enricher(gene = signature_h, universe = background, TERM2GENE = msigdbr_t2g_H)
b <- barplot(ego_h[["SPP1"]], title = "MSigDB H enrichments in SPP1", showCategory = 3 )
tmp8 <- b$data
tmp8$ID <- str_trunc(tmp8$ID,150, "right")
tmp8$ID <- factor(tmp8$ID, levels = tmp8$ID[order(tmp8$Count)])
tmp8$type <- "SPP1"
##CBioPlanet
ego_h[["SPP1"]] <- enricher(gene = signature_h, universe = background, TERM2GENE = CBioPlanet)
b <- barplot(ego_h[["SPP1"]], title = "CBioPlanet enrichments in SPP1", showCategory = 4 )
tmp9 <- b$data
tmp9$ID <- str_trunc(tmp9$ID, 150, "right")
tmp9$ID <- factor(tmp9$ID, levels = tmp9$ID[order(tmp9$Count)])
tmp9$type <- "SPP1"
p7 <- ggplot(tmp7, aes(Count, ID)) +
geom_bar(stat = "identity", color="black", fill = colors_anno["SPP1"]) +
theme_classic()+
geom_text(
aes(label = (paste( "padj=", format(p.adjust,scientific = TRUE, digits = 2)))),
color = "black",
size = 6,
hjust=1,
position = position_dodge(0.5) ) +
theme(text = element_text(size=18)) +
ggtitle("SPP1") +
ylab("")
p8 <- ggplot(tmp8, aes(Count, ID)) +
geom_bar(stat = "identity", color="black", fill = colors_anno["SPP1"]) +
theme_classic()+
geom_text(
aes(label = (paste( "padj=", format(p.adjust,scientific = TRUE, digits = 2)))),
color = "black",
size = 6,
hjust=1,
position = position_dodge(0.5) ) +
theme(text = element_text(size=18)) +
ggtitle("") +
ylab("")
p9 <- ggplot(tmp9, aes(Count, ID)) +
geom_bar(stat = "identity", color="black", fill = colors_anno["SPP1"]) +
theme_classic()+
geom_text(
aes(label = (paste( "padj=", format(p.adjust,scientific = TRUE, digits = 2)))),
color = "black",
size = 6,
hjust=1,
position = position_dodge(0.5) ) +
theme(text = element_text(size=18)) +
ggtitle("") +
ylab("")
p3 <- ggplot(tmp9, aes(Count, ID)) +
geom_bar(stat = "identity", color="white", fill = "white") +
theme_minimal()+
geom_text(
aes(label = (paste( "padj=", format(p.adjust,scientific = TRUE, digits = 2)))),
color = "white",
size = 6,
hjust=1,
position = position_dodge(0.5) ) +
theme(text = element_text(size=18), axis.text = element_text(color = "white")) +
ggtitle("") + theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())+
ylab("BioPlanet") + xlab ("")
(p1 + p4 + p7) + plot_layout(ncol=3)

DT::datatable(rbind(p1$data, p4$data, p7$data),
caption = ("Extended Data Figure 8E"),
extensions = 'Buttons',
options = list( dom = 'Bfrtip',
buttons = c( 'csv', 'excel')))
(p2 + p5 + p8) + plot_layout(ncol=3)

DT::datatable(rbind(p2$data, p5$data, p8$data),
caption = ("Extended Data Figure 8F"),
extensions = 'Buttons',
options = list( dom = 'Bfrtip',
buttons = c( 'csv', 'excel')))
(p3 + p6 + p9) + plot_layout(ncol=3)

DT::datatable(rbind(p3$data, p6$data, p9$data),
caption = ("Extended Data Figure 8G"),
extensions = 'Buttons',
options = list( dom = 'Bfrtip',
buttons = c( 'csv', 'excel')))
Extended Data Figure 8g
IFN_TAM <- list(c("CASP1","CASP4", "CCL2","CCL3","CCL4","CCL7","CCL8", "CD274","CD40","CXCL2","CXCL3","CXCL9","CXCL10","CXCL11", "IDO1","IFI6","IFIT1","IFIT2","IFIT3","IFITM1","IFITM3","IRF1","IRF7","ISG15","LAMP3","PDCD1LG2","TNFSF10","C1QA" ,"C1QC","CD38","IL4I1","ISG15","TNFSF10","IFI44L"))
macro <- AddModuleScore(macro,features = IFN_TAM, name = "IFN_TAM", nbin = 5 )
p <- SCpubr::do_FeaturePlot(macro, features = "IFN_TAM1", pt.size = 0.2, legend.length = 8)
DT::datatable(p$data,
caption = ("Extended Data Figure 8g"),
extensions = 'Buttons',
options = list( dom = 'Bfrtip',
buttons = c( 'csv', 'excel')))
Inflam_TAM <- list(c("CCL2","CCL3","CCL4","CCL5","CCL20","CCL3L1","CCL3L3","CCL4L2","CCL4L4","CXCL1","CXCL2","CXCL3","CXCL5","CXCL8","G0S2","IL1B","IL1RN","IL6","INHBA","KLF2", "KLF6","NEDD9","PMAIP1","S100A8","S100A9","SPP1"))
macro <- AddModuleScore(macro,features = Inflam_TAM, name = "Inflam_TAM", nbin = 5 )
p <- SCpubr::do_FeaturePlot(macro, features = "Inflam_TAM1", pt.size = 0.2, legend.length = 8)
DT::datatable(p$data,
caption = ("Extended Data Figure 8g"),
extensions = 'Buttons',
options = list( dom = 'Bfrtip',
buttons = c( 'csv', 'excel')))
LA_TAM <- list(c("ACP5","AOPE","APOC1","ATF1","C1QA", "C1QB", "C1QC","CCL18","CD163","CD36","CD63","CHI3L1","CTSB","CTSD", "CTSL","F13A1","FABP5","FOLR2","GPNMB","IRF3","LGALS3","LIPA","LPL","MACRO","MerTK","MMP7","MMP9", "MMP12","MRC1","NR1H3","NRF1","NUPR1","PLA2G7","RNASE1","SPARC","SPP1","TFDP2","TREM2","ZEB1"))
macro <- AddModuleScore(macro,features = LA_TAM, name = "LA_TAM", nbin = 5 )
p <- SCpubr::do_FeaturePlot(macro, features = "LA_TAM1", pt.size = 0.2, legend.length = 8)
DT::datatable(p$data,
caption = ("Extended Data Figure 8g"),
extensions = 'Buttons',
options = list( dom = 'Bfrtip',
buttons = c( 'csv', 'excel')))
Angio_TAM <- list(c("ADAM8","AREG","BNIP3","CCL2","CCL4","CCL20","CD163","CD300E","CD44","CD55","CEBPB","CLEC5A","CTSB","EREG","FCN1","FLT1","FN1","HES1","IL1B","IL1RN","CXCL8","MAF","MIF","NR1H3","OLR1","PPARG","S100A8","S100A9","S100A12","SERPINB2","SLC2A1","SPIC","SPP1","THBS1","TIMP1","VCAN","VEGFA"))
macro <- AddModuleScore(macro,features = Angio_TAM, name = "Angio_TAM", nbin = 5 )
p <- SCpubr::do_FeaturePlot(macro, features = "Angio_TAM1", pt.size = 0.2, legend.length = 8)
DT::datatable(p$data,
caption = ("Extended Data Figure 8g"),
extensions = 'Buttons',
options = list( dom = 'Bfrtip',
buttons = c( 'csv', 'excel')))
Reg_TAM <- list(c("CCL2","CD274","CD40","CD80","CD86","CHIT1","CX3CR1","HLA-A","HLA-C","HLA-DQA1","HLA-DQB1","HLA-DRA","HLA-DB1", "HLA-DB5","ICOSLG","IL10","ITGA4","LGALS9","MACRO","MRC1","TGFB2"))
macro <- AddModuleScore(macro,features = Reg_TAM, name = "Reg_TAM", nbin = 5 )
p <- SCpubr::do_FeaturePlot(macro, features = "Reg_TAM1", pt.size = 0.2, legend.length = 8)
DT::datatable(p$data,
caption = ("Extended Data Figure 8g"),
extensions = 'Buttons',
options = list( dom = 'Bfrtip',
buttons = c( 'csv', 'excel')))
Prolif_TAM <- list(c("CCNA2","CDC45","CDK1","H2AFC","HIST1H4C","HMGB1","HMGN2","MKI67","RRM2","STMN1","TOP2A","TUBA1B","TUBB","TYMS"))
macro <- AddModuleScore(macro,features = Prolif_TAM, name = "Prolif_TAM", nbin = 5 )
p <- SCpubr::do_FeaturePlot(macro, features = "Prolif_TAM1", pt.size = 0.2, legend.length = 8)
DT::datatable(p$data,
caption = ("Extended Data Figure 8g"),
extensions = 'Buttons',
options = list( dom = 'Bfrtip',
buttons = c( 'csv', 'excel')))
RTM_TAM <- list(c("C1QA","C1QB","C1QC","CCL7","CD163","CD5L","CD74","CETP","FOLR2","HLA-DPA1", "HLA-DPB1","HLA-DRB1","MACRO","MAF","MS4A7","SLC40A1","VCAM1","VSIG4"))
macro <- AddModuleScore(macro,features = RTM_TAM, name = "RTM_TAM", nbin = 5 )
p <- SCpubr::do_FeaturePlot(macro, features = "RTM_TAM1", pt.size = 0.2, legend.length = 8)
DT::datatable(p$data,
caption = ("Extended Data Figure 8g"),
extensions = 'Buttons',
options = list( dom = 'Bfrtip',
buttons = c( 'csv', 'excel')))
MT_RTM_TAM <- list(c("MT1E","MT1F","MT1G","MT1H","MT1M","MT1X", "MT2A","MTG1E","MTG1M","CXCL9", "CXCL10","GBP1","IFITM3","LGALS2","PARP14","SLC30A1","SLC39A8","STAT1","TNFSF13B"))
macro <- AddModuleScore(macro,features = MT_RTM_TAM, name = "MT_RTM_TAM", nbin = 5 )
p <- SCpubr::do_FeaturePlot(macro, features = "MT_RTM_TAM1", pt.size = 0.2, legend.length = 8)
DT::datatable(p$data,
caption = ("Extended Data Figure 8g"),
extensions = 'Buttons',
options = list( dom = 'Bfrtip',
buttons = c( 'csv', 'excel')))
Classical_TIM <- list(c("S100A4","S100A6","S100A8","S100A9","S100A12","S100A13","AREG","CCL4","CCR2","CD14","CD300E","CD36","CEBPD","CLEC11A","CLEC12A","CSF1R","CTSS","CX3CR1","CXCL2","CXCR4","EGR1","FCN1","FOS","FTH1","IL1B","ITGAM","JUNB","LGALS2","LY6C","LYZ","MAF","MAFB","MPEG1","NFKB1","NFKBIA","NLRP3","NR4A1","NR4A2","OSM","PTGS2","RGS2","SELL","THBS1","VCAN"))
macro <- AddModuleScore(macro,features = Classical_TIM, name = "Classical_TIM", nbin = 5 )
p <- SCpubr::do_FeaturePlot(macro, features = "Classical_TIM1", pt.size = 0.2, legend.length = 8)
DT::datatable(p$data,
caption = ("Extended Data Figure 8g"),
extensions = 'Buttons',
options = list( dom = 'Bfrtip',
buttons = c( 'csv', 'excel')))
NonClassical_Monocytes <- list(c("CX3CR1","CCL5","CD52","CDH23","CDKN1C","FCGR3A", "FCGR3B","FCGR4","GNLY","GZMB","ICAM2","IFITM1","IL32","KLRC4","KLRK1","LILRA", "LILRB1","MS4A7","MTSS1","PRF1","SLC44A2"))
macro <- AddModuleScore(macro,features = NonClassical_Monocytes, name = "NonClassical_Monocytes", nbin = 5 )
p <- SCpubr::do_FeaturePlot(macro, features = "NonClassical_Monocytes1", pt.size = 0.2, legend.length = 8)
DT::datatable(p$data,
caption = ("Extended Data Figure 8g"),
extensions = 'Buttons',
options = list( dom = 'Bfrtip',
buttons = c( 'csv', 'excel')))
Intermediate_Monocytes <- list(c("CCL2","CCL8","CX3CR1","CXCL10","FCGR3A","IFITM2","IL1RN","IRF7","ISG15","ISG20","LILRB1","LILRA1","MS4A7","MTSS1","RHOC","SERPINA1","SIGLEC10","STAT1","TCF7L2","TNFSF10","TNFSF13B"))
macro <- AddModuleScore(macro,features = Intermediate_Monocytes, name = "Intermediate_Monocytes", nbin = 5 )
p <- SCpubr::do_FeaturePlot(macro, features = "Intermediate_Monocytes1", pt.size = 0.2, legend.length = 8)
DT::datatable(p$data,
caption = ("Extended Data Figure 8g"),
extensions = 'Buttons',
options = list( dom = 'Bfrtip',
buttons = c( 'csv', 'excel')))
v1 <- Seurat::VlnPlot(macro,
features = "IFN_TAM1",
group.by = "anno_l1",
pt.size = 0 ,
cols = colors_anno,
) & NoLegend()
v2 <- Seurat::VlnPlot(macro,
features = "Inflam_TAM1",
group.by = "anno_l1",
pt.size = 0 ,
cols = colors_anno,
) & NoLegend()
v3 <- Seurat::VlnPlot(macro,
features = "LA_TAM1",
group.by = "anno_l1",
pt.size = 0 ,
cols = colors_anno,
) & NoLegend()
v4 <- Seurat::VlnPlot(macro,
features = "Angio_TAM1",
group.by = "anno_l1",
pt.size = 0 ,
cols = colors_anno,
) & NoLegend()
v5 <- Seurat::VlnPlot(macro,
features = "Reg_TAM1",
group.by = "anno_l1",
pt.size = 0 ,
cols = colors_anno,
) & NoLegend()
v6 <- Seurat::VlnPlot(macro,
features = "Prolif_TAM1",
group.by = "anno_l1",
pt.size = 0 ,
cols = colors_anno,
) & NoLegend()
v7 <- Seurat::VlnPlot(macro,
features = "RTM_TAM1",
group.by = "anno_l1",
pt.size = 0 ,
cols = colors_anno,
) & NoLegend()
v8 <- Seurat::VlnPlot(macro,
features = "MT_RTM_TAM1",
group.by = "anno_l1",
pt.size = 0 ,
cols = colors_anno,
) & NoLegend()
v9 <- Seurat::VlnPlot(macro,
features = "Classical_TIM1",
group.by = "anno_l1",
pt.size = 0 ,
cols = colors_anno,
) & NoLegend()
v10 <- Seurat::VlnPlot(macro,
features = "NonClassical_Monocytes1",
group.by = "anno_l1",
pt.size = 0 ,
cols = colors_anno,
) & NoLegend()
v11 <- Seurat::VlnPlot(macro,
features = "Intermediate_Monocytes1",
group.by = "anno_l1",
pt.size = 0 ,
cols = colors_anno,
) & NoLegend()
v1 + v2 +v3+v4+v5+v6+v7+v8+v9+v10+v11

SCpubr::do_FeaturePlot(macro, features = "FCGR3A", legend.title = "FCGR3A")
